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The Josephson effect in cychcally coupled Bose-Einstein condensates is studied theoretically. 
We analyze the simultaneous Gross-Pitaevskii equations with coupling terms between adjacent 
condensates. Depending on the initial relative phases between condensates, Josephson current 
flows cyclically to create a quantized vortex. Reducing the coupling between condensates changes 
the motion from periodic to chaotic, thus suppressing the cyclic current. The relation to the 
Kibble-Zurek mechanism is discussed; the density of the generated vortices is less than that 
predicted by this mechanism. 
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The experimental observation of the Bose-Einstein 
condensation (BEC) in a dilute gas of neutral, trapped 
atoms opens up studies of interesting quantum phenom- 
ena on the macroscopic scale .EP The broken gauge sym- 
metry due to BEC produces a macroscopic wave function 
(the order parameter) \l/(r, i) = y/pe^^ , with p as the con- 
densate density and 6 as its phase. This wave function 
obeys the Gross-Pitaevskii equation (GPE), which has 
successfully explained many phenomena of the BEC of 
neutral atoms. One of the significant effects caused by 
this wave function is the Josephson effect. Raghavan et 
al. theoretically studied the coherfint atomic tunneling 
between two weakly coupled BECs.a) They analyzed the 
Feynman two-state model derived from the GPE to in- 
vestigate the characteristic nonlinear phenomena which 
are not accessible with superconductor Josephson junc- 
tions. 

The first motive of this work is to improve the density 
and phase dynamics for multicoupled BECs; this theoret- 
ical study can lead to new experiments. The second mo- 
tive is closely related to the experimental observation of 
vortex generation in superfluid ^He and ''He. Symmetry- 
breaking phase transitions in the early universe imme- 
diately after the Big Bang are expected to leave be- 
hind long-lived topologically stable structures such as 
monopoles, strings or domain walls. Zurek discussed 
the analogy of the second-order phase transition between 
cosmological strings and quantized vortices in the su- 
perfluid, and suggested a cryogenic experiment which 
tested the cosmological scenario for string formation.H-' 
The Kibble-Zurek (KZ) scenario is the following: after 
the system enters the broken symmetry phase through 
a rapid phase transition, the ordered phases start to 
form, simultaneously and independently, in many parts 
of the system. Each ordered region has an independent 
phase 9 of the order parameter field ^'(r,i), because 
they are disconnected causally from each other. Critical 
slowing down near the transion temperature suppresses 
the velocity with which the order parameter coherence 



spreads. Consequently, when the disconnected domains 
grow together and superfluid coherence becomes estab- 
lished throughout the entire volume, the topological de- 
fects around which supercurrent flows, i.e., the vortices, 
are frozen out at the boundaries. Zurek estimated the 
defect density from the coherent length which is com- 
parable to the frozen domain size, it should be noted 
that this estimation is based on the assumption that the 
boundary between more than three domains always cre- 
ates a defect. This theoretical scenario iiujtivated various 
experimental works. The experimentaj'EF which heated 
small regions in superfluid ^He by a neutron-induced nu- 
clear reaction yielded results consistent with the theoret- 
ical estimation of the vortex tangle density, while Volovik 
subseqently noted that the experimental workfl) using ro- 
tating superfluid '^He evidenced not only the vortices gen- 
erated by the KZ_mechanism but also those induced by 
the rotating flowfl^ The experimentQ* which expanded 
liquid ■*He rapidly through the lambda line observed a 
density that was at least two orders of magnitude less 
than the theoretical estimation; the reason for the dis- 
agreement remains unresolved. Hence, at present, the 
KZ scenario is not necessarily confirmed experimentally. 

Thus, this paper discusses one important factor which 
may be omitted in the above scenario. A cyclic super- 
current should flow through adjacent weakly coupled do- 
mains in order for a vortex to appear there. However, 
it is not clear whether cyclically coupled domains nec- 
essarily result in cyclic Josephson current (CJC) or not. 
Therefore, by improving the physics of two weakly cou- 
pled BECs, this work studies three cyclically coupled 
BECs. Whether CJC appears or not is found to depend 
on the initial relative phases and the coupling constants 
between adjacent BECs. It is impossible to control the 
initial relative phases. Hence, if initial relative phases 
are assumed to take arbitrary values between and 2t: 
with equal probability, we obtain the expected value of 
CJC which depends on the coupling constants. By ap- 
plying these results to more than three cyclically coupled 
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BECs, wc find that the density of generated vortices may 
be reduced compared with the original KZ scenario. 

First, the formulation_for the two weakly coupled 
BECs is described briefly .B The macroscopic wave func- 
tion \E'(r,t) in a trap potential V{r) at T — satisfies 
the GPE 



Ft 



-[V{r)+go\^{r,t)\^]^{r,t), 



where the cyclic boundary condition means "00 = V'« 
and ipn+i = "01 ■ The parameter Kij+i is the coupling 
constant between ith and i-l- 1th BEC. The simplest case 
of n = 3 is discussed. The insertion of V'i = V^s^^ in 
eq. (7) yields 

nNi = 2Ki2y/NiN2sm(j)i2 - 2/^31 V'TVgTV^ sin ,^3(8a) 

hN2 = 2/^23 y^Va^Va sin 023 - 2is:i2 v/^V^sin0i^8b) 



(1) 

with go = 47r?i^a/m, m the atomic mass, and a the s- 
wave scattering length of the atoms. When V{r) is a 
double-well trap, ^'(r,t) is written as 



2if3iV^V3]Vi sin03i - 2/4:23 V^VaiVg sin 02:(8c) 
-El - UiNi 



+K12 J cos 012 + i^3i 1 



(2) 



where $1^2 describe the condensate which is localized in 
each trap, "01,2(0 = ■\/^i,2e*^^'^^*-' and the total num- 



ber of atoms Nt — Ni + N2 — \ipi 



1^-2 



The weak 




coupling causes the overlap integral / ^i{r)^2{''')dr to 
be negligible, which leads to the nonlinear two-mode dy- 
namical equations 



+K23\I ^ cos 4123 + Ki2^ 
I\2 



COS 031 -I- K; 





■ COS 031 



■ COS 012 



COS 023, 



(8d) 



dt 

.9l/>2 



[El + UiNi)^i - K^2 



iTl-^ = {E^2 + C^2A^2)^2 - i^^i 

at 



(3a) 
(3b) 



Ml 

M2 = 
Ms = 

+1 

where 0ij = 9i — 6j is a relative phase. The Josephson 
supercurrent which flows from the ith to the jth conden- 
sate is expressed as 

— ^y]V^sin0y 



(8f) 



J, 



(9) 



Equations (8a)- (8f) conserve the total number iV^ 
A^i + A^2 + ^3 and the total energy 



Here, Ei 2 are the zero-point energies in each well, 
C^i,2-^i,2 the atomic self-interaction energies and K the 
amplitude of the tunneling between condensates. In- 
troducing the population imbalance z(t) = {Ni{t) — ^ 
N2it))/NT and the relative phase 0(i) = 6'2(t) - 6*1 (i), E = ^{E,N, + -U^Nf)+'^2K,J^yJ^~N~ cos(j),j. (10) 
and rescaling to a dimensionless time 2Kt/h t, 
eqs.(3a) and (3b) are reduced to the canonical form: 

dz dH d(j) dH 



Since the wave function ^'(r,i) is single- valued, the rel- 
ative phases must satisfy 



dt d(t> ' 
Here, the Hamiltonian is 



dt 



dz 



j Xlipdt ■ 



h2 + 023 + 'i>3i = 27rm, 



(11) 



A;: 



H = ^Ez + -^-Vl 



(5) 



and the parameters are 
^ ^ {Ui + U2)Nt El , [Ui - U2)Nt 



4K 



2K 



AK 

(6) 

Although the exact solutions of eq. (Q) are studied in 
detail in Ref.2, investigating the conserved energy H as 
a function of z{t) and 0(t) qualitatively reveals its dy- 
namics. One characteristic motion is "macroscopic quan- 
tum self-trapping" (MQST); the average population dif- 
ference between two BECs is nonzero. Another is tt- 
phase oscillations which are observed experimentally ia 
the weak link separating superfluid '^He-B reservoirs. B 
These Josephson tt states and their quantum decay have 
also been studied by Hatakenaka.E|) 

Extending the above formulation, we will study n 
cyclically coupled BECs. The substitution of the rep- 
resentation ^(r,i) — X^i V'i(0*i'i('") into eq. (1) yields 

dip 

= {E.-L + UiNi)'ilji - ifi^j_i-0j_i - Kij+iip^+i, 

(z = l,2,...,nl7) 



where C is a contour that passes cyclically through the 
condensates and m an integer. The variable 03i is elim- 
inated by eq. (|ll|). Equations (8a)- (8f) show that the 
dynamics is independent of m. The following calcula- 
tions are performed for the dimensionless variables: 

Kt N, hJ,^j . 



Nt 



KNt 



We will study the dynamics for Ei = E2 = E3 = E' 
and Ki — K2 — K3 — K. Then, the parameters for 
the dimensionless equations are e = E' /K , and A^ — 
UiNT/K that represents the strength of the coupling. 
We assume Ai = A2 = A3 = A. Attention will be paid 
mainly to whether CJC flows or not. First, under A = 1 
and the fixed initial values ni(0) = 712(0) = 713(0) = 1/3, 
the initial phases 0i2(O) and 023(0) are changed. The 
oscillation of ji^j is shown in Fig. 1. In Fig. 1(a) for 
012(0) — 023(0) = 7r/8, each current oscillates about 
zero; net currents flow neither clockwise nor counter- 
clockwise. However, when the initial phases are changed 
(Figs. 1(b) and 1(c)), the time-averaged value of every 
current becomes negative, which means that CJC ap- 
pears. 
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(a) A=l : <A,/O)=;r/8,0^^(O)=;r/8 



(b) A=l : «!',,(0)=2;r/3 , ,*^^(0)=7;r/l2 
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Fig. 1. Oscillation of ji^2 (solid line), J2^3 (dashed-dotted line) 
and J3_,i (dotted line). Values of A and initial relative phases 
are shown in each figure. 



The appearance of CJC is found by investigating the 
time- and piece-averaged current 



< J >< 



1 r 



+ J2^3+J3^l), (13) 



where T is a period of time which is sufficiently longer 
than the period of oscillation. 

Figure 2(a) shows how < j >av depends on 0i2(O) 
and (/'23(0) at A = f . The current < j >av has two ex- 
tremums for ^12(0) — 023(0) — 27r/3 and 47r/3 which 
give the stationary solutions. However, when A is in- 
creased, i.e., the coupling is reduced, as shown in Figs. 
1(b) and 1(d), the oscillation of currents become chaotic; 
the transition from periodic to chaotic motion is con- 
firmed by the analysis of the power spectrum. The 
important property of the n > 3 BECs dynamics is 
this transition to chaos which is absent in two coupled 
BECs.B When the motion becomes chaotic, < j >av is 
lessened considerably. Thus, it is concluded that the ap- 
pearance of CJC depends on the initial relative phases 
and the coupling constants, which are summarized in 
Fig. 2. Since we cannot control the initial relative 
phases, CJC appears with some probability that depends 
on the coupling constants. If every relative phase is 
assumed to take a value between and 27r with equal 
probability, the expected values or the ensemble aver- 
ages of < j >av for A =1, 3 and 10 are, respectively, 
0.174, 0.088 and 0.029. This system has some unique 
modes which are compared with those of two coupled 
BECs.B' The oscillative ones are zero-phase, 27r/3-phase 
and 47r/3-phase modes in which the time-average of the 
phase < (pij >— 0, 27r/3 and 47r/3, respectively. There 
are also running-phase modes that create MQST. The 
analysis of the linear stability finds that the zero-phase 
mode is stable for an arbitrary value of A, while the 27r/3- 
and 47r/3-phase modes become unstable for A > 9/4. 





Fig. 2. Dependence of the time-averaged current < j >av on the 
initial relative phases 012(0) and </!>23(0) for (a) A =1, (b) 3 and 
(c) 10. 



grow up after the quench of the system, as described in 
the introduction of this paper, the coupling constants 
between adjacent domains are expected not to be con- 
stant but to increase with time. Such time-dependent 
coupling constants change the motion gradually from ini- 
tially chaotic to periodic. In a real system, every pair of 
condensates need not have the same coupling constant 



When the independent domains of the ordered phase ^^ven some dissipative term is introduced to eq. 
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(1), the calculation for K12 > K23 = K31 shows that the 
first and second condensates are united, subsequently the 
united condensate and the third one create the two-body 
motion; the dissipative term is described in the following 
paper. 

The above consideration is applicable to n cyclically 
coupled BECs. This system has 27rfc/7i-phase modes 
about the stationary solutions which result from the ini- 
tial relative phases ^12(0) = (j)23{Q) = • • • = (/)„_i,„(0) = 
27rfc/n, where fc = 0, 1, • • • , n — 1. In order to obtain the 
exact conclusions on CJC in n BECs needs to carry out 
an analysis for each value of n. However, the following 
estimation for weak coupling roughly gives the depen- 
dence of CJC on n. Following the case of n = 3, we will 
define the ensemble-averaged CJC for n BECs 

< JCJC >n= ^^J^n-l j\<J (^)M^ (14) 

X = (012(0), (/.23(0), • • • , K-l,nm. (15) 

where 

< j >a« (33) ^ ^ jfc^fc+i (i, a;) j dt, (16) 

is the time- and piece-averaged current starting from the 
initial phase x. When the coupling is weak, as shown in 
Fig. 2(c), the main contribution for CJC comes from the 
stationary solutions. In the limit of weak coupling, we 
assume that only the stationary solutions produce CJC: 



< jcjc > 




It is easy to show that < jcjc >n decreases with 
n to disappear for n ^ 00. For example, we find 

< JCJC >n= 2.9 X 10-2,4.0 X 10-3,7.8 x IQ-* from eq. 
(17), respectively, for n = 3,4, 5. Even if the coupling is 
not weak, the relative ratio < jcjc >n / < jcJC >3 is 
not changed very much. When independent domains are 
nucleated in a system, the number of n cyclically cou- 
pled BECs decreases with n. Hence it is concluded that 
CJC arises mainly from cyclic BECs with a small value 
of n. Following these results, we can estimate how the 
vortex density in the KZ mechanism is reduced. Then, 
the first assumption is that the vortices are generated 
only from the CJC of n = 3 cyclic BECs based on 
the above considerations. The second is that the vor- 
tex density in the KZ mechanism is caused by the CJC 
due to the stationary solution 0i2 ~ 023 — 27r/3, i.e., 

< JKZ >3— 0.577 which is the maximum of < jcjc >3j 
indicating that CJC flows independently of the initial rel- 
ative phases. However, the CJC actually depends on the 
initial relative phases and the coupling constant, so that 

< jcJC >3 is less than < jxz >3 to lessen the density 
of generated vortices. Figure 3 shows the dependence of 

< jcJC >3 / < JKZ >3 on A, which can be related to 
the reduced ratio of the vortex density. 
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Fig. 3. Dependence of < jcjC >3 / < jKZ >3 on A. 



studies that include the application to a real system and 
dissipative effects are to be published shortly. 

The authors thank T. lida for useful discussions. 
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It is difhcult to quantitatively estimate the coupling 
constants in superfluid helium, although it may be pos- 
sible in a dilute gas of trapped atoms. More detailed 



